Objetivo General

Description: This code aims to download the FGJ Victims Database using datos.cdmx.gob.mx, Mexico City’s open data portal. Using districts provided by IECM: https://datos.cdmx.gob.mx/dataset/coloniascdmx

Setup .-

dev.off()
rm(list=ls())
# install.packages(pacman) #Se utiliza para manejar paquetes de R
pacman::p_load(st,sf,rgeos,tidyverse,jsonlite,rgdal,leaflet,ggspatial,raster,htmltools)

Mapa IECM - 2019 .-

Vamos a usar el mapa del IECM que se encuentra en https://datos.cdmx.gob.mx/dataset/coloniascdmx

getwd()
## [1] "/Users/dortega/Desktop/ADIP/social tic/output"
mapa_iecd <- readOGR("/Users/dortega/Desktop/ADIP/social tic//input/mgpc_2019.shp") #Cargamos el mapa de colonias del IEECM
## OGR data source with driver: ESRI Shapefile 
## Source: "/Users/dortega/Desktop/ADIP/social tic/input/mgpc_2019.shp", layer: "mgpc_2019"
## with 1815 features
## It has 8 fields
class(mapa_iecd)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"
mapa_iecd@data %>% head()
##   ENT CVEDT        NOMDT DTTOLOC  CVEUT                  NOMUT POB2010 ID
## 0   9     2 AZCAPOTZALCO      05 02-001               AGUILERA    2014  1
## 1   9     2 AZCAPOTZALCO      05 02-002                 ALDANA    3378  2
## 2   9     2 AZCAPOTZALCO      05 02-005          ANGEL ZIMBRON    2737  3
## 3   9     2 AZCAPOTZALCO      05 02-006                 ARENAL    4817  4
## 4   9     2 AZCAPOTZALCO      03 02-007 CENTRO DE AZCAPOTZALCO    3043  5
## 5   9     2 AZCAPOTZALCO      05 02-008               CLAVERIA    9570  6
#Ver el Mapa
plot(mapa_iecd, col="#f2f2f2", bg="skyblue", border="black" ,lwd=0.5 )

Víctimas en Carpetas de Investigación de la CDMX .-

Consulta base de datos de víctimas de la fiscalía con la API; podemos encontrar aquí la documentacion de la API de CKAN: https://docs.ckan.org/en/2.9/api/

Pero la construcción de una consulta requiere esencialmente de:

  • URL de consulta
  • función que se pide
  • id de Base de Datos que necesitamos

En nuestro caso también añadiremos dos parámetros de consulta adicionales:

  • Límite
  • Offset
url <- "https://datos.cdmx.gob.mx/api/3/action/"
funcion <- "datastore_search?resource_id="
id_carpetas <- "d543a7b1-f8cb-439f-8a5c-e56c5479eeb5" #carpetas y viictimas de investigación en la CDMX:
limite <- "&limit=32000" # Si no se especifica el número de observaciones se regresan 100 por default
# limita las consultas de json a 32,000 observaciones por consulta la página de datos abiertos por el momento

Hacemos el request, para ello hacemos uso de fromJSON, ya que la API de CKAN contesta con archivos JSON, esta función le da estructura de lista para su manejo.

req <- fromJSON(paste0(url,funcion,id_carpetas,limite))
class(req) # el request es una lista, por lo que debemos buscar lo que pedimos dentro de la lista
## [1] "list"
glimpse(req) # objetos de la lista con listas y objetos dentro de ellos
## List of 3
##  $ help   : chr "https://datos.cdmx.gob.mx/api/3/action/help_show?name=datastore_search"
##  $ success: logi TRUE
##  $ result :List of 10
##   ..$ include_total             : logi TRUE
##   ..$ limit                     : int 32000
##   ..$ records_format            : chr "objects"
##   ..$ resource_id               : chr "d543a7b1-f8cb-439f-8a5c-e56c5479eeb5"
##   ..$ total_estimation_threshold: NULL
##   ..$ records                   :'data.frame':   32000 obs. of  23 variables:
##   .. ..$ _id            : int [1:32000] 1 2 3 4 5 6 7 8 17 9 ...
##   .. ..$ idCarpeta      : int [1:32000] 8324429 8324430 8324431 8324435 8324438 8324442 8324444 8324451 8324482 8324454 ...
##   .. ..$ Ano_inicio     : int [1:32000] 2019 2019 2019 2019 2019 2019 2019 2019 2019 2019 ...
##   .. ..$ Mes_inicio     : chr [1:32000] "Enero" "Enero" "Enero" "Enero" ...
##   .. ..$ FechaInicio    : chr [1:32000] "2019-04-01T00:00:00" "2019-04-01T00:00:00" "2019-04-01T00:00:00" "2019-04-01T00:00:00" ...
##   .. ..$ Delito         : chr [1:32000] "FRAUDE" "PRODUCCIÓN, IMPRESIÓN, ENAJENACIÓN, DISTRIBUCIÓN, ALTERACIÓN O FALSIFICACIÓN DE TÍTULOS AL PORTADOR, DOCUMENTOS"| __truncated__ "ROBO A TRANSEUNTE SALIENDO DEL BANCO CON VIOLENCIA" "ROBO DE VEHICULO DE SERVICIO PARTICULAR SIN VIOLENCIA" ...
##   .. ..$ Categoria      : chr [1:32000] "DELITO DE BAJO IMPACTO" "DELITO DE BAJO IMPACTO" "ROBO A CUENTAHABIENTE SALIENDO DEL CAJERO CON VIOLENCIA" "ROBO DE VEHÍCULO CON Y SIN VIOLENCIA" ...
##   .. ..$ Sexo           : chr [1:32000] "Masculino" "Femenino" "Masculino" "Masculino" ...
##   .. ..$ Edad           : chr [1:32000] "62" "38" "42" "35" ...
##   .. ..$ TipoPersona    : chr [1:32000] "FISICA" "FISICA" "FISICA" "FISICA" ...
##   .. ..$ CalidadJuridica: chr [1:32000] "OFENDIDO" "VICTIMA Y DENUNCIANTE" "VICTIMA Y DENUNCIANTE" "VICTIMA Y DENUNCIANTE" ...
##   .. ..$ competencia    : chr [1:32000] "FUERO COMUN" "FUERO COMUN" "FUERO COMUN" "FUERO COMUN" ...
##   .. ..$ Ano_hecho      : chr [1:32000] "2018" "2018" "2018" "2019" ...
##   .. ..$ Mes_hecho      : chr [1:32000] "Agosto" "Diciembre" "Diciembre" "Enero" ...
##   .. ..$ FechaHecho     : chr [1:32000] "29/08/2018" "15/12/2018" "22/12/2018" "04/01/2019" ...
##   .. ..$ HoraHecho      : chr [1:32000] "12:00:00" "15:00:00" "15:30:00" "06:00:00" ...
##   .. ..$ HoraInicio     : chr [1:32000] "12:19:00" "12:20:00" "12:23:00" "12:27:00" ...
##   .. ..$ AlcaldiaHechos : chr [1:32000] "ALVARO OBREGON" "AZCAPOTZALCO" "COYOACAN" "IZTACALCO" ...
##   .. ..$ ColoniaHechos  : chr [1:32000] "GUADALUPE INN" "VICTORIA DE LAS DEMOCRACIAS" "COPILCO UNIVERSIDAD ISSSTE" "AGRÍCOLA PANTITLAN" ...
##   .. ..$ Calle_hechos   : chr [1:32000] "INSUGENTES SUR" "AV.  CUATLAHUAC" "COPILCO" "CALLE 6" ...
##   .. ..$ Calle_hechos2  : chr [1:32000] "NA" "NA" "NA" "ENTRE PRIVADA DEL VALLE Y PRIVADA GONZALEZ" ...
##   .. ..$ latitud        : chr [1:32000] "19.36125" "19.47181" "19.33797" "19.40327" ...
##   .. ..$ longitud       : chr [1:32000] "-99.18314" "-99.16458" "-99.18611" "-99.05983" ...
##   ..$ fields                    :'data.frame':   23 obs. of  2 variables:
##   .. ..$ id  : chr [1:23] "_id" "idCarpeta" "Ano_inicio" "Mes_inicio" ...
##   .. ..$ type: chr [1:23] "int" "numeric" "numeric" "text" ...
##   ..$ _links                    :List of 2
##   .. ..$ start: chr "/api/3/action/datastore_search?resource_id=d543a7b1-f8cb-439f-8a5c-e56c5479eeb5&limit=32000"
##   .. ..$ next : chr "/api/3/action/datastore_search?resource_id=d543a7b1-f8cb-439f-8a5c-e56c5479eeb5&limit=32000&offset=32000"
##   ..$ total                     : int 679500
##   ..$ total_was_estimated       : logi FALSE
req$success # el request fue exitoso
## [1] TRUE

Habiendo revisado que la consulta fue exitosa, extraemos la información del requiest enviado:

#Extraemos los datos de la consulta
datos_fgj <- req$result$records
glimpse(datos_fgj) #Revisamos qué datos vienen
## Rows: 32,000
## Columns: 23
## $ `_id`           <int> 1, 2, 3, 4, 5, 6, 7, 8, 17, 9, 10, 11, 12, 13, 14, 15,…
## $ idCarpeta       <int> 8324429, 8324430, 8324431, 8324435, 8324438, 8324442, …
## $ Ano_inicio      <int> 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, …
## $ Mes_inicio      <chr> "Enero", "Enero", "Enero", "Enero", "Enero", "Enero", …
## $ FechaInicio     <chr> "2019-04-01T00:00:00", "2019-04-01T00:00:00", "2019-04…
## $ Delito          <chr> "FRAUDE", "PRODUCCIÓN, IMPRESIÓN, ENAJENACIÓN, DISTRIB…
## $ Categoria       <chr> "DELITO DE BAJO IMPACTO", "DELITO DE BAJO IMPACTO", "R…
## $ Sexo            <chr> "Masculino", "Femenino", "Masculino", "Masculino", "Ma…
## $ Edad            <chr> "62", "38", "42", "35", "NA", "42", "55", "13", "NA", …
## $ TipoPersona     <chr> "FISICA", "FISICA", "FISICA", "FISICA", "FISICA", "FIS…
## $ CalidadJuridica <chr> "OFENDIDO", "VICTIMA Y DENUNCIANTE", "VICTIMA Y DENUNC…
## $ competencia     <chr> "FUERO COMUN", "FUERO COMUN", "FUERO COMUN", "FUERO CO…
## $ Ano_hecho       <chr> "2018", "2018", "2018", "2019", "2019", "2018", "2019"…
## $ Mes_hecho       <chr> "Agosto", "Diciembre", "Diciembre", "Enero", "Enero", …
## $ FechaHecho      <chr> "29/08/2018", "15/12/2018", "22/12/2018", "04/01/2019"…
## $ HoraHecho       <chr> "12:00:00", "15:00:00", "15:30:00", "06:00:00", "20:00…
## $ HoraInicio      <chr> "12:19:00", "12:20:00", "12:23:00", "12:27:00", "12:35…
## $ AlcaldiaHechos  <chr> "ALVARO OBREGON", "AZCAPOTZALCO", "COYOACAN", "IZTACAL…
## $ ColoniaHechos   <chr> "GUADALUPE INN", "VICTORIA DE LAS DEMOCRACIAS", "COPIL…
## $ Calle_hechos    <chr> "INSUGENTES SUR", "AV.  CUATLAHUAC", "COPILCO", "CALLE…
## $ Calle_hechos2   <chr> "NA", "NA", "NA", "ENTRE PRIVADA DEL VALLE Y PRIVADA G…
## $ latitud         <chr> "19.36125", "19.47181", "19.33797", "19.40327", "19.35…
## $ longitud        <chr> "-99.18314", "-99.16458", "-99.18611", "-99.05983", "-…
total_records <- req$result$total #Vemos cuál es el total de los datos
total_records
## [1] 679500

Como podemos apreciar, el número de datos de la base es menor al extraído, por lo que necesiutamos extraer la base secuencialmente hasta terminar de extraer el archivo, lo haremos con el siguiete loop:

while (nrow(datos_fgj)<total_records){
  ofst <- nrow(datos_fgj)
  ofset <- paste0("&offset=",ofst)
  req_2 <- fromJSON(paste0(url,funcion,id_carpetas,limite,ofset))
  datos_fgj <- rbind(datos_fgj,req_2$result$records)
  print(nrow(datos_fgj))
}
## [1] 64000
## [1] 96000
## [1] 128000
## [1] 160000
## [1] 192000
## [1] 224000
## [1] 256000
## [1] 288000
## [1] 320000
## [1] 352000
## [1] 384000
## [1] 416000
## [1] 448000
## [1] 480000
## [1] 512000
## [1] 544000
## [1] 576000
## [1] 608000
## [1] 640000
## [1] 672000
## [1] 679500
#limpiamos espacio
rm(url,funcion,id_carpetas,limite,ofset,ofst,req,req_2,total_records)
glimpse(datos_fgj)
## Rows: 679,500
## Columns: 23
## $ `_id`           <int> 1, 2, 3, 4, 5, 6, 7, 8, 17, 9, 10, 11, 12, 13, 14, 15,…
## $ idCarpeta       <int> 8324429, 8324430, 8324431, 8324435, 8324438, 8324442, …
## $ Ano_inicio      <int> 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, …
## $ Mes_inicio      <chr> "Enero", "Enero", "Enero", "Enero", "Enero", "Enero", …
## $ FechaInicio     <chr> "2019-04-01T00:00:00", "2019-04-01T00:00:00", "2019-04…
## $ Delito          <chr> "FRAUDE", "PRODUCCIÓN, IMPRESIÓN, ENAJENACIÓN, DISTRIB…
## $ Categoria       <chr> "DELITO DE BAJO IMPACTO", "DELITO DE BAJO IMPACTO", "R…
## $ Sexo            <chr> "Masculino", "Femenino", "Masculino", "Masculino", "Ma…
## $ Edad            <chr> "62", "38", "42", "35", "NA", "42", "55", "13", "NA", …
## $ TipoPersona     <chr> "FISICA", "FISICA", "FISICA", "FISICA", "FISICA", "FIS…
## $ CalidadJuridica <chr> "OFENDIDO", "VICTIMA Y DENUNCIANTE", "VICTIMA Y DENUNC…
## $ competencia     <chr> "FUERO COMUN", "FUERO COMUN", "FUERO COMUN", "FUERO CO…
## $ Ano_hecho       <chr> "2018", "2018", "2018", "2019", "2019", "2018", "2019"…
## $ Mes_hecho       <chr> "Agosto", "Diciembre", "Diciembre", "Enero", "Enero", …
## $ FechaHecho      <chr> "29/08/2018", "15/12/2018", "22/12/2018", "04/01/2019"…
## $ HoraHecho       <chr> "12:00:00", "15:00:00", "15:30:00", "06:00:00", "20:00…
## $ HoraInicio      <chr> "12:19:00", "12:20:00", "12:23:00", "12:27:00", "12:35…
## $ AlcaldiaHechos  <chr> "ALVARO OBREGON", "AZCAPOTZALCO", "COYOACAN", "IZTACAL…
## $ ColoniaHechos   <chr> "GUADALUPE INN", "VICTORIA DE LAS DEMOCRACIAS", "COPIL…
## $ Calle_hechos    <chr> "INSUGENTES SUR", "AV.  CUATLAHUAC", "COPILCO", "CALLE…
## $ Calle_hechos2   <chr> "NA", "NA", "NA", "ENTRE PRIVADA DEL VALLE Y PRIVADA G…
## $ latitud         <chr> "19.36125", "19.47181", "19.33797", "19.40327", "19.35…
## $ longitud        <chr> "-99.18314", "-99.16458", "-99.18611", "-99.05983", "-…

Ahora que tenemos la base de datos completa podemos elegir las variables que necesitamos y operar sobre ellas.

Datos de Interés.-

Para el tutorial de hoy queremos únicamente la categoría del delito, el año del hecho y el sexo de la victima, aunque las demás variables pueden ser de amplio interés para otros análisis.

datos_fgj <- datos_fgj %>% tibble() %>%  dplyr::select(Categoria,Sexo,Ano_hecho,latitud,longitud) # elegimos lo que queremos 
glimpse(datos_fgj)
## Rows: 679,500
## Columns: 5
## $ Categoria <chr> "DELITO DE BAJO IMPACTO", "DELITO DE BAJO IMPACTO", "ROBO A …
## $ Sexo      <chr> "Masculino", "Femenino", "Masculino", "Masculino", "Masculin…
## $ Ano_hecho <chr> "2018", "2018", "2018", "2019", "2019", "2018", "2019", "201…
## $ latitud   <chr> "19.36125", "19.47181", "19.33797", "19.40327", "19.3548", "…
## $ longitud  <chr> "-99.18314", "-99.16458", "-99.18611", "-99.05983", "-99.063…
#Problema, necesitamos que las latitudes y longitudes sean numéricos:
datos_fgj$latitud <- as.numeric(datos_fgj$latitud) #hago números las coordenadas
## Warning: NAs introducidos por coerción
datos_fgj$longitud <- as.numeric(datos_fgj$longitud) #números de las coordenadas
## Warning: NAs introducidos por coerción
sum(is.na(datos_fgj$latitud))
## [1] 28542
sum(is.na(datos_fgj$longitud))
## [1] 28546
datos_fgj <- datos_fgj[-which(is.na(datos_fgj$latitud)),] # quitamos los que no tienen latitud
datos_fgj <- datos_fgj[-which(is.na(datos_fgj$longitud)),] # quitamos los que no tienen longitud

Vemos los puntos de latitud y longitud que extraímos:

ggplot(data = datos_fgj,aes(longitud,latitud)) +
  geom_point(size=0.01) # datos normales

Para los subsecuentes análisis es importante ver la frecuencia temporal de las observaciones:

table(datos_fgj$Ano_hecho)
## 
##   1950   1952   1955   1958   1962   1963   1967   1968   1969   1971   1972 
##      1      2      2      3      4      3      2      5     18      2      8 
##   1973   1974   1975   1976   1979   1980   1981   1982   1983   1984   1985 
##      2      4      1      5      2      5      8      6     10     12      4 
##   1986   1987   1988   1989   1990   1991   1992   1993   1994   1995   1996 
##      8      6      5      9      5      7      8      5     13     19     17 
##   1997   1998   1999   2000   2001   2002   2003   2004   2005   2006   2007 
##      7     32     20     49     54     59     26     45     43     59     87 
##   2008   2009   2010   2011   2012   2013   2014   2015   2016   2017   2018 
##    128    151    263    271    334    354    638    774   1317   2919  20127 
##   2019   2020   2021     NA 
## 287013 186153 149651    169
plot(table(datos_fgj$Ano_hecho))

Vamos a restringir datos a los hechos entre 2019 - 2021, dado que es la fecha donde comienzan a capturarse estos datos en este formato:

datos_fgj <- datos_fgj %>% filter(Ano_hecho > 2018) %>% filter(Ano_hecho!= "NA")
plot(table(datos_fgj$Ano_hecho))

Gráficas Estáticas.-

Densidad de crimen por concentración de puntos:

ggplot(datos_fgj, aes(x = longitud, y = latitud)) + 
  coord_equal() + 
  xlab('Longitud') + 
  ylab('Latitud') + 
  labs(title="Densidad de crimen en la CDMX con datos entre 2019-2021")+
  stat_density2d(aes(fill = ..level..), alpha = .5,
                 h = .02, n = 300,
                 geom = "polygon", data = datos_fgj) + 
  scale_fill_viridis_c() 

Víctimas en carpetas por sexo:

ggplot(data = datos_fgj,aes(longitud,latitud)) +
  geom_point(size=0.01,aes(colour=Sexo))+
  labs(title="Carpetas de Investigación por Sexo de la Víctima en la CDMX entre 2016 - 2021",
       caption = "Con información de la FGJ en Datos Abiertos: https://datos.cdmx.gob.mx/dataset/victimas-en-carpetas-de-investigacion-fgj")

Quitamos NA:

datos_sexo <- datos_fgj[!datos_fgj$Sexo=="NA",]
nrow(datos_fgj)-nrow(datos_sexo) # Datos perdidos al quitar NA
## [1] 130173
ggplot(data = datos_sexo,aes(longitud,latitud)) +
  geom_point(size=0.01,aes(colour=Sexo)) + 
  theme(legend.background = element_rect(fill="lightblue",size=0.5, linetype="solid", colour ="darkblue")) +
  labs(title="Carpetas de Investigación por Sexo de la Víctima en la CDMX entre 2016 - 2021",
       caption = "Con información de la FGJ en Datos Abiertos: https://datos.cdmx.gob.mx/dataset/victimas-en-carpetas-de-investigacion-fgj")

Víctimas porCategoría del delito:

ggplot(data = datos_fgj,aes(longitud,latitud)) +
  geom_point(size=0.01,aes(colour=Categoria)) + 
  theme(legend.position = "bottom",
        legend.background = element_rect(fill="lightblue",size=0.5, linetype="solid", colour ="darkblue")) +
  guides(colour=guide_legend(ncol=3))+
  labs(title="Carpetas de Investigación por Categoría de Delito en la CDMX entre 2016 - 2021",
       caption = "Con información de la FGJ en Datos Abiertos: https://datos.cdmx.gob.mx/dataset/victimas-en-carpetas-de-investigacion-fgj")

Dificil de interpretar, podemos seccionarlo, un mapa por cada categoría del delito —

ggplot(data = datos_fgj,aes(longitud,latitud)) +
  geom_point(size=0.01,aes(colour=Categoria)) + 
  facet_wrap( ~ Categoria)+
  theme(legend.position="none")+
  labs(title="Carpetas de Investigación por Categoría de Delito en la CDMX entre 2016 - 2021",
       caption = "Con información de la FGJ en Datos Abiertos: https://datos.cdmx.gob.mx/dataset/victimas-en-carpetas-de-investigacion-fgj")

Intersectar Mapas-

Clase de ambos objetos es distinta:

class(mapa_iecd)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"

Necesitamos hacer que los datos de latitud y longitud sean coordinadas:

#Extraemos los puntos y hacemos una matriz de ellos:
puntos_fgj <- matrix(c(datos_fgj$longitud,datos_fgj$latitud),ncol = 2)
puntos_fgj %>% head()
##           [,1]     [,2]
## [1,] -99.18314 19.36125
## [2,] -99.16458 19.47181
## [3,] -99.18611 19.33797
## [4,] -99.05983 19.40327
## [5,] -99.06324 19.35480
## [6,] -99.16016 19.33537

Ambos deben tener la misma predicción al convertir en objeto espacial Si revisamos el diccionario de datos encontrramos la projección utilizada por la FGJ: https://datos.cdmx.gob.mx/dataset/victimas-en-carpetas-de-investigacion-fgj
Podemos obtener el prj4 de internet: https://spatialreference.org/ref/epsg/wgs-84/

#forzamos a hacer espaciales
puntos_fgj <- SpatialPoints(puntos_fgj,proj4string = CRS("+proj=longlat +datum=WGS84 +no_defs"),bbox=NULL) 
#proj4string se utiliza para obtener la proyección de un objeto espacial
proj4string(mapa_iecd)
## Warning in proj4string(mapa_iecd): CRS object has comment, which is lost in output; in tests, see
## https://cran.r-project.org/web/packages/sp/vignettes/CRS_warnings.html
## [1] "+proj=utm +zone=14 +datum=WGS84 +units=m +no_defs"
#Traducimos al CRS (Coordinate Reference System) del mapa
puntos_fgj <- spTransform(puntos_fgj,crs(mapa_iecd))
puntos_fgj
## class       : SpatialPoints 
## features    : 650954 
## extent      : 464033.8, 505585.6, 2114926, 2165379  (xmin, xmax, ymin, ymax)
## crs         : +proj=utm +zone=14 +datum=WGS84 +units=m +no_defs
mapa_iecd
## class       : SpatialPolygonsDataFrame 
## features    : 1815 
## extent      : 463188.4, 505712.8, 2114935, 2166535  (xmin, xmax, ymin, ymax)
## crs         : +proj=utm +zone=14 +datum=WGS84 +units=m +no_defs 
## variables   : 8
## names       : ENT, CVEDT,          NOMDT, DTTOLOC,  CVEUT,       NOMUT, POB2010,   ID 
## min values  :   9,    10, ALVARO OBREGON,      01, 02-001, 10 DE ABRIL,    1000,    1 
## max values  :   9,     9,     XOCHIMILCO,      33, 17-086,  ZOTOLTITLA,    9998, 1815
# Ahora ambos son objetos espaciales en las mismas coordenadas y projecciones
class(puntos_fgj) 
## [1] "SpatialPoints"
## attr(,"package")
## [1] "sp"
class(mapa_iecd)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"
st_crs(puntos_fgj)==st_crs(mapa_iecd)
## [1] TRUE

Ahora veamos el mapa superpuesto:

plot(mapa_iecd, bg="lightblue", main="Delitos en la CDMX entre 2016-2020")
points(puntos_fgj,pch=20,col="blue")

La diferencia principal entre los mapas de antes y ahora, es que ahora ambos objetos son espaciales: por lo que ahora podemos operar sobre ellos. Hacemos uso de la función over para revisar cuántos puntos hay sobre cada polígono.

over(mapa_iecd,puntos_fgj,returnList = FALSE, fn=sum) %>% head()
##    0    1    2    3    4    5 
## 2203 7613 2668  469  145  407
#guardamos la variable
mapa_iecd@data$delitos_totales <- over(mapa_iecd,puntos_fgj,returnList = FALSE, fn=sum)

Leaflet.-

Para poder visualizar estos datos de una manera agradable tornamos a leaflet; sin embargo debemos de cambiar de nuevo la proyección a la que este último software utiliza de manera default:

class(mapa_iecd)
## [1] "SpatialPolygonsDataFrame"
## attr(,"package")
## [1] "sp"
mapa_iecd <- spTransform(mapa_iecd, "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
## Warning in Polygon(coords = crds): less than 4 coordinates in polygon

Hacemos las divisiones con los deciles de delitos totales y con los mismos el color del mapa:

bins <- quantile(mapa_iecd@data$delitos_totales,probs = seq(.1, 1, by = .1),na.rm=T)
bins <- as.numeric(levels(factor(bins)))
pal <- colorBin("YlOrRd", domain = mapa_iecd$delitos_totales, bins = bins)
m <- leaflet(mapa_iecd) %>%
  addTiles()

m %>% addPolygons(data=mapa_iecd, weight = .5,  fillColor = ~pal(delitos_totales), 
              opacity = 1,
              color = "black",
              fillOpacity = 0.7,
              label = labels,
              labelOptions = labelOptions(
                style = list("font-weight" = "normal", padding = "3px 8px"),
                textsize = "15px",
                direction = "auto"))

labels <- sprintf(
  "<strong>%s</strong><br/>%g Delitos totales entre 2019-2021",
  mapa_iecd$NOMUT, mapa_iecd$delitos_totales) %>% lapply(htmltools::HTML)

m